3D Medical Imaging Page 1 THREE DIMENSIONAL MEDICAL IMAGING Michael W. Vannier Daniel Geist Donald E. Gayou Mallinckrodt Institute of Radiology Washington University School of Medicine 510 S. Kingshighway Blvd. St. Louis, Mo. 63110 ABSTRACT Surface and volumetric three dimensional imaging methods have found application in diagnostic medical imaging and research. The acquisition, modeling, classification, and computer graphics rendering of discrete image volumes will be introduced. Applications in diagnosis (craniofacial, orthopedic, cardiovascular, and others) as well as reconstruction methods for generic serial sections will be described. C-Software for 3D reconstruction which operates on an IBM PC/AT clone is described. INTRODUCTION The substantial quantity of information generated by CT and MRI scanners has presented both a problem and an opportunity for clinicians. Vast quantities of data must be processed before diagnostically relevant information can be extracted. On the other hand, steady progress in the development of data reduction and visualization software has created a class of clinician-users who desire to manipulate imaging data to provide views of the data compatible with pertinent diagnostic and therapeutic needs in the management of a particular patient. For example, a radiologist may wish to apply various image enhancement techniques to better visualize small details a difficult case while a surgeon may wish to interactively manipulate the 3-D data on screen for surgical planning. The radiation oncologist may utilize CT or MRI image data in a slightly different way to facilitate true 3-D radiation therapy planning. In each of these cases, physicians with differing requirements wish to manipulate the same 3-D data. In such situations, a workstation architecture offers an ideal solution to the problem of distributed medical imaging. Networking permits multiple users to efficiently access imaging data stored at a central location. The price/performance ratio of workstations is steadily improving, while cost per seat continues to decline. Increasingly powerful chipsets utilizing ASIC (Application-specific integrated circuits) provide powerful graphics engines for complex image processing and graphical manipulation of large quantities of 3-D data. Recently, the distinction between engineering workstations and high end PC's has blurred, and medical imaging software packages have become available for PC's at reasonable cost. There are a number of potential advantages inherent in a workstation approach to medical imaging. The user may display and manipulate 3-D data in whatever manner is most appropriate for the task at hand. Direct interaction with the data allows the clinician to comprehend needed information in minimum time. The flexibility and power of this approach is useful not only in clinical settings but also for research and teaching. Finally, trends in both hardware and software technology ensure that the cost of workstations will decline as functionality and "user-friendliness" continue to improve, making these systems available to an ever-broadening group of users. Image processing and computer graphics have become important software tools in the day to day practice of medicine. The use of computers in diagnosis of internal morphological abnormalities is commonplace in radiology departments. X-rays, ultrasound, radioactive drugs - inhaled or injected, and proton magnetic resonance are the physical basis of computerized medical imaging instruments, often called "scanners". Various types of scanners now exist, enabling the extraction of data in the form of images from the human body. Image processing techniques are used to acquire this raw data, reconstruct and enhance it to give the physician maximum information, and as an aid in diagnosis, treatment planning, and evaluation. Three-dimensional (3-D) reconstruction of images has emerged as an important new tool for medical imaging(1). The original data obtained from a scanner is often a set of 2-dimensional contiguous slices (Figure 1) of the body (Figure 2). The reconstruction program takes this data set and renders a 3-D image of internal or external body surfaces. This helps physicians understand the complex 3-D anatomy present in those slices. The objective is to make the reconstruction as useful as possible for the physician. This requires the extraction of as much surface detail as possible in the 3-D image. 3-D reconstructions are used in a variety of medical applications, such as planning of craniofacial surgery(2), neurosurgery(3), and orthopedics(4). Presently, 3-D reconstruction programs run on computed-tomographic (CT) scanner control mini-computers (DEC PDP-11 or Data General Eclipse) or on specialized computer graphics workstations. These units are more costly and are not so widely available as personal computers. We have developed 3-D surface reconstruction software for CT scans in the C language on a DEC Vaxmate personal computer, an IBM PC/AT clone. We believe that an implementation of 3-D reconstruction software on personal computers will popularize the use of these methods for CT scan data, and encourage the use of these reconstructions as a standard tool for physicians who need information on a patient's internal anatomy. Many approaches to the computerized reconstruction of serial slice data sets have been described and implemented(1). These descriptions typically describe the implementation in obscure or incomplete fashion, and lack sufficient detail to allow the user to reproduce the work without an inordinate effort. This paper presents a limited PC implementation with a detailed description of the methods used, and availability of the source code (in the C language). COMPUTED TOMOGRAPHY A CT scan is an X-ray image of the body formed by computerized reconstructions. These slices can vary in thickness from 1 to 15 mm. For high resolution studies the slices are typically 2 mm thick. The pixel values of the images typically represent the average X-ray attenuation of body matter in a corresponding 1 mm x 1 mm x 2 mm cube of the body (which is called a volume element or voxel). This value is obtained by exposing the body to X-ray radiation from several directions while measuring the attenuation of the X-rays passing through the body(4). A projection in one direction is obtained with a moving X-ray source and a set of detectors. By taking all the rays that pass through one voxel in the body and using a method called filtered back-projection(5), one can accurately compute the average x-ray attenuation of that voxel. The x-ray attenuation values computed for each voxel are expressed in a scale relative to water called the Hounsfield scale(6) in honor or G. Hounsfield who (along with A. Cormack) was awarded the Nobel prize in 1979 for inventing the CT scanner. A CT scan image data set consisting of many contiguous slices is a three dimensional array of Hounsfield values (Figure 2). Each CT scan slice is a two dimensional array (Figure 1). The slice we typically use contain 256 x 256 voxels but slice data can be 512 x 512 in size (Figure 2). The number of slices used for one patient's examination typically varies between fifty and one hundred total slices. This results in more than 5 million voxels or 10 MB of data (one voxel is two bytes) for each examination. The manipulation and interpretation of such large amounts of data is a significant problem in medical imaging. THREE DIMENSIONAL RECONSTRUCTION Interpretation of CT slice images requires special training and expertise, and is typically performed by a diagnostic radiologist (M.D. with special training in diagnostic medical imaging). Non-radiologist physicians and surgeons may have difficulty understanding the complex anatomy represented by these two-dimensional (2-D) images. It is often desirable to synthesize a more familiar 3-D image from a sequence of 2-D slices to aid in the interpretation and utilization of CT scan examinations. Methodologies for performing this synthesis have been developed over the past 10 years, and are generally referred to as "3-D reconstruction from serial slices". Similar approaches have been used in electronic microscopy, ultrasound and magnetic resonance imaging. Three-dimensional (3-D) CT scan surface reconstructions were originally developed to simplify the interpretation and improve the utility of CT scans of the face and skull. These CT scans contain morphologic information regarding the skull and soft tissue structures of the head. Unfortunately, the number, complexity and redundancy of CT slice images limits their value in patients with craniofacial disorders. The skull itself is intrinsically three dimensional, and the multiplicity of two-dimensional slice images that result from a CT scan examination is a significant obstacle for the interpreter who must form a complete understanding of complex abnormalities. In the initial development of three-dimensional surface reconstruction methods, the emulation of a dry skull appearance from a set of CT scans was sought. Among physicians and surgeons, the dry skull is a common reference for describing and communicating the results of imaging studies. The similarity between three-dimensional CT scan surface reconstructions and a dry skull is not coincidental. This simulated dry skull appearance for 3-D CT scan reconstructions is a familiar and consistent format that enhances their utility in a clinical setting. These methods have been applied to CT scan surface reconstruction for craniofacial disorders. Three-dimensional surface views are computed from the desired perspective, based on the original sequence of thin CT scan slices. This step typically involves prior knowledge of the orientation for the desired views. Real-time interaction and display have been impractical to date; however, because we use reconstructed images only hours or days after processing, this has not been a significant limitation. The advantages of simplicity and efficiency make processing on modest computer hardware, such as that found in CT scanner control computers, practical. Three dimensional reconstruction from serial CT scan slices can be performed by thresholding the data, projecting the contours and then shading the results. Bone has much higher CT density (x-ray attenuation) than soft tissue or air, therefore the Hounsfield number for voxels containing bone is significantly higher than those of soft tissue. By discarding all image data below some predetermined value, we are left only with voxels which contain bone. In this way we can obtain an image of the bone surface. This technique is called thresholding or level slicing(2). By using an even lower threshold value we can discard only voxels containing air, thus obtaining an image of soft tissue or the outer surface of the patient. Thresholding between various types of soft tissue is not as effective since the CT data does not have such distinct values separating various types of soft tissue. Other types of scanners such as magnetic resonance imaging or MRI scanners allow better differentiation between soft tissues which facilitates 3-D display of their surfaces. After thresholding has been completed, rendering of bone surface is achieved by shading(7). First, a viewing direction is chosen. The 3-D view in that direction is projected onto a 2-D plane called the image plane for eventual display on a CRT screen (Figure 3). Each pixel in the reconstruction image is given a gray scale value according to a suitable mathematical model. The surface points that are visible in the view direction are projected onto the plane of the CRT display screen. There have been many different shading techniques proposed in recent years(8-13). Some of them involve initially transforming the voxel data to a 3-D geometric database(11, 12) such as surface contours or polygons and then creating images from that database. This transformation is an expensive procedure because of the high dimensionality of CT data sets and the complexity of the surfaces they contain. Other methods create the image directly from voxel data(8-10, 13). In this case, we compare two methods, Distance and Gradient shading. DISTANCE SHADING To create distance shades the program passes through the CT data, computing a visible surface image from a specified direction(14). Each voxel is examined individually. If the voxel lies above a given threshold it's distance is saved. The distance is the number of voxels traversed from a fixed viewing site in the chosen direction. Typically it is just the array index (in the specified direction) of the viewer from the voxel where the transition occurred. The other two indices of the voxel are the indices of the corresponding pixel on the image plane. Only the first instance of threshold transition is saved for each pixel in the image plane. In this way we create a Z-buffer(7) of the view in the desired direction. After scaling these distances for some frame buffer and projecting them on a two dimensional display, the image obtained is a view of the object from the chosen direction. The intensity of the projected pixel on the screen is given in the formula: I(i,j)= Imax (D-d)/D Here d is the distance of the voxel from the origin and D is the maximum distance a voxel can be (Figure 3) from the viewing plane. We have included here a small program in C to do depth shading (Distances.c) of an object in the X direction. The CT data is contained in a set of files with the name CTBILD.nnn where nnn is a number. Each file contains one CT slice and the first file is conventionally CTBILD.001 and so on. Distance (or distance-only as it is sometimes called) shading is the fastest method of shading since the computational effort to obtain a surface approximation is minimal. The CT data is read and processed only once and the surface shade is obtained immediately by simply saving the index values of the original slice array when a threshold transition is encountered. However, such shading gives a very smooth appearance that is lacking in fine detail and curvature information. Distance shading was included here primarily for comparison, although in some medical cases it can still give sufficient information. The program (Distances.C) demonstrates the simplest type of depth shaded 3D reconstruction. Today there are systems which can generate much more impressive images. Nevertheless, this small program illustrates the fundamentals of 3-D reconstruction from serial slices. The Distances.C program assumes that the original CT scans are in consecutive files named CTBILD.XXX where XXX are the slice numbers in ascending order. The program reads and processes a slice at a time and produces a line of the image output for each. The result is a shaded projection of the original slice along the viewing direction. The projection is computed by examining each voxel along the viewing direction until the value passes a certain threshold and saving the distance in the corresponding element of the projection vector. Each incoming slice is stored in a matrix in memory. The program codes the array index of the location where the threshold is exceeded, as the grayscale value of the output image. The result gives us the depth-coded 3-D reconstruction image of the scanned object as seen from the chosen view direction. GRADIENT SHADING The reconstruction program takes a set of CT files from one consecutive scan sequence and outputs a set of 3-D view images. The set of output images consists of 9 basic views: Front, back, right-lateral, left-lateral, left midsagittal, right midsagittal, and top intracranial. For each view direction two images are created - one gradient shaded and the other depth coded (or distance only shaded). The views were reconstructed first by thresholding (or level slicing) the data and then shading. A Z-buffer algorithm was used to give the distance-only shaded images (distances.c). The gradient of the bone surface at a point is obtained by examining the distances calculated at the surface for adjacent points along a specified viewing direction. This is a relatively accurate approximation--the error is the size of the unit squared (which in this case is Dx or Dy, the size of a pixel). The surface is asumed to be a continuous, smooth function of (x,y) in a local area around the point. The gradient is computed as ÄS={[Z(x+1,y)-Z(x-1,y)]/2*Dx,[Z(x,y+1)-Z(x,y-1)]/2*Dy,1}. The gradient value computed for each pixel is Cosf where f is the angle between the surface gradient vector and the viewing direction. Gradient shading is exceptionally good in edge rendering of surfaces because at the edges the gradient becomes perpendicular to the viewing direction and the corresponding shade becomes very dark. The rapid changes in gradient intensity at edges amplifies these features in images. This quality of gradient shading enables the depiction of fine anatomical details of internal body surfaces. We interpolated between slices to overcome the problem of artifactual contour creation at slices stacked on one another. The slices were interpolated linearly according to the zoom or magnification factor (the ratio of voxel distance in the z direction and the x direction when the scans are presumed to be parallel to the XY plane). The zoom factor is not necessarily an integer and the number of interpolation points between two slices is not uniform, for all slices pairs. It is important that an image has the correct dimensions to properly scale the results and allow one to make measurements from a life-size hardcopy (photograph) of the final image. Our experiments have shown that doing interpolation on the image gray shades after 3D reconstruction has little effect on image quality. On the other hand, interpolation before the reconstruction gives a substantial improvement in picture quality. The pictures obtained after slice interpolation were improved in apparent surface detail, though they were not completely free of the slice-edge contour artifact. The reconstruction program read the input data twice: forward and backward. The first pass created the rear, left and bottom views and the second pass created the front, right and top views. The method for computing the gradient shade was simple: whenever a threshold transition was encountered, the four columns surrounding the voxel were found. The gradient value was calculated as explained above using the transition direction to be the surface direction (that is, if the transition was in the X direction the surface was assumed to be S=X(y,x)). Flow-charts for the more complex generation of gradient shaded images in am IBM PC/AT clone are provided for a main program which first uses getdistances to obtain depth-coded images that are gradient shaded in the x, y, and z directions using doviews. IMPLEMENTATION An IBM PC/AT clone was used for 3-D reconstruction processing. This unit has a built-in Ethernet connection and software support which enables the transfer of high data volume from VAX computers through a DECNET network. The ability to manipulate high volume of data is the principal problem in medical image 3-D reconstruction processing and display. Typically, the CT scan data set for a single patient will contain 10 MBytes of data. A networked PC/AT is an attractive alternative for CT scan post-processing applications. The PC/AT clone, running on a clock of 8 Mhz, was equipped with one 1Mbyte of memory, a 1.2 MB floppy disk drive and a monochrome display. An Ethernet controller (3 COM) was installed and DECNET local area network (LAN) software was used. This enables the PC/AT clone to transfer files back and forth from a VAX or use the VAX disks as if they were local PC/AT devices (VAX/VMS services for MS-DOS software product from DEC). A 20 MB hard disk was installed, leaving room for two PC/AT bus 16-bit application boards. An additional 2 MB of memory was eventually installed inside the PC/AT clone, as well as an 80287 floating point coprocessor. We did use these options in the latest stages of our development and found that both can be used to decrease the computation time for 3D surface reconstructions. The frame buffer system we used was an AT&T Targa24 board. The Targa images are RGB with 8 bits for each color. Similar results can be obtained with alternative frame buffer interface cards for the PC/AT. The code was written in Microsoft C version 5, but with a few changes such as file handling calls and huge array declarations, it can be compiled by any C compiler. We found the Microsoft CodeView debugger to be an extremely powerful tool in this effort which dramatically reduced our development time. The program is organized in several modules (see flow charts): GETDISTANCES - This is actually the first pass over the CT data. The program travels over the CT data and extracts the distance of the first threshold transition in the X, Y, and Z directions simultaneously. The distance of a pixel on the image plane in the X and Y directions is obtained after processing one slice. While the distance in the Z direction is only known after processing all the slices. The distances are saved in three files respectively as floating point data. (Getdistances flow chart) DOVIEWS - After the distance file exists, this procedure goes over the distance data and creates the distance-only and gradient view. The pixels at the frame of the view are missing some neighbors, so the central difference can not be used to calculate their gradient. Instead forward and backward difference is used. This is why there is a separate code to handle the first and last line of the view, and the first and last pixel of each line. (Doviews flow chart) The patient used in the example reconstruction is a three an a half year old girl with a cleft lip and palate. A total of 78 CT scan slices were used in the reconstruction, a total of 5,111,808 voxels. The program which produced the 3-D views ran an hour and twenty minutes total or about one minute for each slice. If enough RAM memory is available on the system, a RAM-disk can be defined and all the intermediate distance values can be saved in memory instead of on the hard disk. This method can dramatically reduce processing time for surface reconstruction. CONCLUSION We have demonstrated the feasibility of 3-D image reconstruction on a personal computer, and our work shows that satisfactory images can be obtained in reasonable time. The newest generation of personal computers, such as the 80386 or 68020 baed machines, can be used as the basis of a medical image workstation and perform many useful post-processing imaging functions. We have found that a workstation will need at least 5 MB of memory and a 40 MB hard disk to enable convenient utilization. However, for those who have an IBM PC/AT or compatible, good quality 3-D images can be obtained with modest effort. Computer graphics have an important role in biomedical research and medical imaging. The emergence of high performance graphics workstations and medical image transmission display and archiving systems are important new developments. 3-D reconstruction from CT and MRI scans to meet the clinical requirements for evaluating patients with complex morphologic abnormalities can be accomplished using available methods. These methods have been applied to diagnostic and therapeutic problems in orthopedics, cardiac and craniofacial surgery, and nuclear medicine. The basis for analysis of complex craniofacial malformations has been static surface views of the skull and facial structures in predefined projections. More elaborate analysis, using a CAD/CAM system, has been applied to especially complex cases. REFERENCES Geist D, Vannier MW, Knapp RH: Networked personal computer software system for 3-D CT scan reconstruction imaging, National Computer Graphics Assn., Vol. III, 113-117, 1988. Geist D, Vannier MW, Knapp RH, Marsh JL: 3-D medical imaging in an object oriented environment, SPIE Medical Imaging II, Vol. 914, 1299-1306, 1988. Herman GT, Liu HK: Three-dimensional display of human organs from computed tomographs, Comp. Graphics and Image Proc., 7:130-138, 1978. Herman GT, Liu HK: Three dimensional display of human organs from computed tomographs, Jour. Computer Assisted Tomography, 1:155-160, 1977. Hildebolt CF, Vannier MW: 3-D measurement accuracy of skull surface landmarks, American Journal of Physical Anthropology, 76:497-503, 1988. Hildebolt CF, Vannier MW, Knapp RH: Validation of skull 3-D CT measurements, American Journal of Physical Anthropology, in press. Requicha AAG: Representations for rigid solids: Theory, methods, and systems, Comp. Surveys, 12:(4)437-464, Dec. 1980. Srihari SN: Representation of three-dimensional digital images, Comp. Surveys, 13:(4)399-424, 1981. Totty WG, Vannier MW: Analysis of complex musculoskeletal anatomy using three-dimensional surface reconstruction, Radiology, 150:173-177, 1984. Udupa JK: Display of 3D information in discrete 3D scenes produced by computerized tomography, Proc. IEEE, 71:(3)420-431, 1983. Vannier MW, Marsh JL, Warren JO: Three-dimensional computer graphics for craniofacial surgical planning and evaluation, Comput. Graphics, 17:263-274, 1983. Vannier MW, Marsh JL, Gado MH: Three-dimensional display of intracranial soft-tissue abnormalities, American Journal of Neuroradiology, 4:520-521, 1983. Vannier MW, Gutierrez FR, Laschinger J, Gronemeyer S, Canter C, Knapp RH: Three-dimensional MR imaging of congenital heart disease, Radiographics, Vol. 8, No. 5, 857-872, 1988. Zatz LM: Basic principles in computed tomography scanning, in T.H. Newton and D.G. Potts, Eds: Radiology of the Skull and Brain-Technical Aspects of Computed Tomography, 5:3853-3876, (1981) C.V. Mosby, St. Louis.